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Resonant graphene dopants, such as hydrogen adatoms, experience long-range effective interaction 
mediated by conduction electrons. As a result of this interaction, when several adatoms are present 
in the sample, hopping of adatoms between sites belonging to different sublattices involves significant 
energy changes. Different inelastic mechanisms facilitating such hopping - coupling to phonons and 
conduction electrons - are considered. It is estimated that the diffusion of hydrogen adatoms is 
rather slow, amounting to roughly one hop to a nearest neighbor per millisecond. 

PACS numbers: 


I. INTRODUCTION 


Graphene is a two-dimensional materiaP^ holding a 
lot of technological potential. Diffusion of hydrogen on 
graphene is a problem of major interest for a number of 
applications. One potential application arises from the 
possibility of using graphene for hydrogen storage. For its 
successful implementation it is important to understand 
and be able to predict the rates of hydrogen absorption, 
desorption, and diffusion on graphene. 

A different class of applications concerns graphene elec¬ 
tronics. One of the main obstacles here is graphene’s 
good electric conduction in the intrinsic state. It is 
thus desirable to develop ways of suppressing graphene 
conductivity and turning graphene into a semiconduc¬ 
tor in a controllable way. Avenues explored to achieve 
this objective include a number of possibilities: open¬ 
ing a gap in graphene bilayers with an interlayer 
bia^^HD^ applying elastic strairP®^, carving out finite- 
width nanoribbon^l^Ifl^, inducing strong spin-orbital 
couplin^isHIll^ or using chemical dopin^i^HIll 

Hydrogen is one especially promising dopant. Gom- 
plete coverage of graphene with hydrogen atoms, how¬ 
ever, results in a dielectric (graphane) with a very large 
gap, ~ 5eV, see Refs. [231 HOI a situation equally unfavor¬ 
able for electronics applications. Nonetheless, partially 
hydrogenated graphene remains a viable candidate. But 
for an incomplete coverage of graphene with hydrogen, 
the question of diffusion becomes important. This ques¬ 
tion is made much more interesting and non-trivial by 
the existence of an effective interaction between the hy¬ 
drogen atoms. Let us briefly explain the origin of this 
interaction before discussing how it might affect hydro¬ 
gen transport. 

A hydrogen impurity is resonant; its spectrum has an 
energy level close to the Dirac point of the conduction 
TT-band of graphen^^. Resonant hopping of conduction 
electrons on and off the hydrogen level gives rise to a large 
scattering amplitude of conduction electrons. It turns out 
that resonant scattering leads to an effect quite similar to 
that of a vacancy (a lattice site rendered unaccessible to 
TT-electrons because of the lack of a carbon atom there). 


or to a very strong potential of a substitution defect. As 
a result, the wave functions of electrons are significantly 
modified. This leads to a long-range interaction between 
dopants mediated by conduction electrons. In its origin, 
this interaction is similar to the RKKY interaction or 
the classic Gasimir effect mediated by virtual photonP^. 
One notable difference with the RKKY interaction is that 
the latter is usually considered for distances exceeding 
the Fermi wavelength, kpR ^ 1, while in the case of 
intrinsic graphene one has to deal with the opposite limit, 
kpR ^ 1. 

For vacancies and also for strong substitution defects 
or resonant impurities, the interaction energy W (i?) of 
two impurities can already be estimated from dimen¬ 
sional considerations. Graphene 7r-electrons have a gap¬ 
less Dirac spectrum, E = ±vp. Thus, for infinitely 
strong impurities in an intrinsic graphene, one can con¬ 
struct only a single combination with the dimension of 
energy: W{R) ~ hv/R. This interaction was first found 
in Ref. |26| for kp = 0. Its dependence on the chemical 
potential p = vpp was elucidated in Ref. H?) Due to 
the presence of two sublattices in a honeycomb graphene 
lattice and the resulting quantum interference, the sign 
of the interaction depends on whether the dopants reside 
on the same (HA-case) or opposite sublattices (HR-case). 
In the same-sublattice case the interaction is repulsive. 


Waa^ 


Unv cos^Oaa 
2R\n^{R/a)' 


( 1 ) 


The phase angle 9aa (R-) = cos cj) depends on both 
the length of the inter-impurity radius-vector R and the 
angle (j) it makes with a zigzag direction, see Fig. 1. 

The interaction between two impurities residing on dif¬ 
ferent sublattices is more interesting since it is sensitive 
to the chemical potential p: 

Wab{R) = - 2 fR 0 (fR - \t\), ( 2 ) 

2itln {R/a) 

where 6^45 (R) = cos </>-I- </>. The first term in Eq. (j 2 | 
is repulsive. Its origin is the same as that of Waa- the 
band spectrum of the conduction electrons is modified by 
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FIG. 1: Graphene inter-atomic distance is a = 0.14nm. 
Atoms belonging to different sublattices (A and B) are repre¬ 
sented with different colors (red and blue). Two on-site impu¬ 
rities (grey circles) are placed on graphene (AB-configuration 
shown). Periodic boundary conditions are assumed in the 
armchair (x) and zigzag (y) directions. The angle (f) is counted 
from a zigzag direction (the t/-axis). 


the interaction with the two impurities to an extent that 
depends on their relative positions. The second term is 
the result of the formation of two bound states, one above 
and one below the Dirac point, at E = ±i5r, with the 
absolute value of the energy given bjEl 

c _ 7iu|sin6(4B| 

^ R\n{R/a)' 

The second term in Eq. § simply represents the en¬ 
ergy of an occupied bound state, —2 £r,, with the factor 
2 accounting for the spin degeneracy. This term is log¬ 
arithmically dominant (for In {R/a) ^ 1) over the first 
term in Eq. but only when the chemical potential 
is confined between the upper and lower bound states, 
^ ^ ^R,. However, when the chemical potential 

moves to either below or above the two levels, both of 
them become empty or filled, respectively, and their con¬ 
tribution to Wab disappears. This is accounted for by 
the 0-function in Eq. ([^. The sign of the interaction of 
the two dopants can thus be changed by a mere variation 
of the chemical potential^. 

This pairwise interaction of hydrogen atoms, Eqs. Q- 
([^ , may have a profound effect on the migration tenden¬ 
cies in favor of either clustering or spreading, depending 
on the chemical potential, which controls the attractive 
or repulsive nature of the interaction. It is thus impor¬ 
tant to evaluate the migration rate of hydrogen adatoms. 

Let us start with the propagation of a single hydro¬ 
gen atom on an ideal graphene crystal. Hydrogen atoms 
reside above carbon atoms. Because of the large mass 
of the hydrogen atom (compared with the electron mass 
for example) the tunneling rate between carbon sites is 
rather small but not insignificant, and we are going to 
show that it would still result in a rapid spreading of the 
adatoms. 

However, we are also going to see that the presence of 
a second adatom, even many interatomic distances away. 


changes the situation dramatically. Indeed, the differ¬ 
ence of the energies before and after tunneling, AE, is of 
the order of several to tens of meV. This energy is many 
orders of magnitude larger than the broadening of the 
hydrogen level due to elastic tunneling. Therefore, in or¬ 
der to occur, a tunneling event must be assisted by some 
mechanism susceptible to deliver or carry away the en¬ 
ergy difference. The two candidates for this are phonons 
and electrons, so we are going to calculate the rate of 
such phonon-assisted and electron-assisted tunneling of 
hydrogen on graphene. 

Throughout the paper we utilize the effective Dirac 
model of graphene spectrum applicable at low energies. 
Correspondingly, our results are valid when typical dis¬ 
tances between adatoms are much larger than the carbon 
lattice spacing. For such long distances the numerical 
methods, such as DFT, would be impractical. On the 
other hand, our approach cannot be applied for short dis¬ 
tances, in particular to the problem of hy drogen dimers 
where numerical methods are needecpS129| 

The plan of the paper is as follows. In Section [H] 
we review the hopping of a single hydrogen adatom on 
graphene. We find the inelastic tunneling rates for the 
processes involving phonons and electrons, in Sections 
HI and W respectively. Section |V] is a discussion of the 


results. 


II. HYDROGEN HOPPING ON GRAPHENE 

First-principle calculations of the height of the po¬ 
tential barrier between first neighbors for hydrogen 
adatoms resulted in a broad range of valuea^SHSSl from 
0.29eVtol.3eV. This is a consequence of the relatively 
small size of the system simulated by DFT combined with 
the important role played by the carbon lattice relaxation 
in the binding of adatoms. 

Coronene and coronene-like species have been shown to 
effectively simulate the electronic structure of graphene 
hydrogenation, especially near central carbon sites®. Us¬ 
ing a development version of the Q-Chem 4.3 quantum 
chemistry packag^®, we modeled the potential landscape 
of a hydrogen adatom above a central carbon atom of cir- 
cumcoronene -\- H). All DFT calculations were 

carried out on the doublet electronic ground state us¬ 
ing the B3LYP functional with the double-zeta polarized 
basis set 6-31G(d). 

A full geometry optimization with the adatom was first 
computed allowing carbon lattice relaxation. After lat¬ 
tice relaxation, DFT results show an in-plane carbon- 
carbon distance a = 0.14 nm with a torsion angle of 
9 = 15 degrees. The hydrogen adatom was then scanned 
along an evenly-spaced square 0.7 nm x 0.7nm grid above 
the carbon atom allowing relaxation of the z-component 
of the gradient for the hydrogen (see Fig. and con¬ 
straining the carbon lattice to that of the optimized 
structure. The minimum potential energy was taken at 
each step. 
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FIG. 2: Hydrogenated circumcoronene {C^ 4 ,H\^s + H), the 
adsorption site is shown for reference. 


In order to get an approximation for the hopping am¬ 
plitude, we truncated the potential obtained from DFT 
to within an in-plane disk of radius a/2, keeping the po¬ 
tential constant outside of this region. This allows us to 
set up the following Hamiltonian in the position repre¬ 
sentation 

where Mh is the mass of a hydrogen atom and 
UoFTi^’ y) truncated DFT potential. The ground 

state wave function of the Hamiltonian Q extended to 
large distances r a/2 was obtained using the Fourier 
grid Hamiltonian methocP^. This gives us an approxima¬ 
tion to the evanescent tail of the wave function in that 
region. Due to the large size of the grid, diagonalization 
was done numerically using the Primme iterative diago¬ 
nalization routin^^. The results are shown in Fig. 

Using this model of a ground state wave function, 
we numerically estimate the hopping amplitude t = 
—0.61 /reV by integrating the product of two neighbor¬ 
ing well ground state wave functions multiplied by the 
DFT potential over a circular region of radius a/2 cen¬ 
tered on one of the two sites. At the same time, we can 
estimate the overlap integral for the wave functions of 
hydrogen states residing on adjacent atoms: /r-r' = 
/ dr'0*(r — R)'i/i(r — R'), which plays a determinant role 
in the evaluation of the transition rate. Using the same 
model for the ground state wave function as above, nu¬ 
merical integration gives /r_ = 3.8 X 10 This ap¬ 

proximation of /r-r' is an underestimate by construc¬ 
tion of the model ground state wave function, as it ig¬ 
nores the other wells in the lattice. However, the small 
values of both t and /r-r' are reflecting the sharply lo¬ 
calized nature of the hydrogen wave functions illustrated 
in Fig.[^ a posteriori justifying the approximation made. 


A single hydrogen atom hopping over an ideal graphene 
sheet must therefore behave like a band particle and 
propagate with the band velocity u ~ at/Ti ^ 1 cm/s. 
This velocity is much smaller than the band velocity of 
TT-electrons, u = 1 x 10® m/s. Nonetheless, if hydro¬ 
gen atoms absorbed by graphene were moving with such 
velocity u, this would result in a rather quick homoge¬ 
nization of their distribution. 

However, such a small (by microscopic scales) velocity 
means that the hydrogen bandwidth is very narrow. This 
should cause the suppression of elastic hopping as soon as 
the on-site potential energies between neighboring sites 
are different by more than t. This is why the interaction 
mediated by conduction electrons makes elastic hopping 
impossible. When one of the atoms hops, it virtually 
always changes the sublattice, A —>■ H or H —)■ A: the 
likelihood of tunneling to a second-nearest neighbor on 
the same sublattice site is negligible, because of the much 
larger distance that must be covered. Suppose now that 
a second hydrogen atom happens to be some distance 
R away, sitting atop a carbon atom belonging to the 
same sublattice. The interaction between the atoms then 
results in a positive energy Waa(^) > 0. When one 
of the atoms hops to its nearest neighbor carbon atom, 
it lands on the opposite sublattice, and the interaction 
energy in the final state is negative, IU 4 B(R+a) < 0. The 
total energy change in this process, Waa{^) — Wab(^ + 
a), therefore, exceeds both Waa and Wab- Even for two 
atoms sitting as far as a thousand interatomic distances 
apart, the energy change is estimated from Eqs. Q-© 
to be of the order of 1 meV. This energy is many orders 
of magnitude greater than the hydrogen adatom hopping 
bandwidth. We conclude, therefore, that the interaction 
mediated by conduction electrons makes elastic hopping 
impossible and results in a collective pinning of hydrogen 
on graphene. 

Consequently, the tunneling can occur only if facil¬ 
itated by additional processes which supply (or carry 
away) the energy difference. Such processes must in¬ 
volve excitations of graphene, phonons or electron-hole 
pairs. The transition rates in both cases are evaluated in 
Sections ITH] and |IV| respectively. 


III. HYDROGEN-PHONON INTERACTION 

Let us consider inelastic tunneling of a hydrogen atom 
between adjacent sites R and R' where the atom has 
potential energy Hr and Wr', respectively. The extra 
energy Hr — Hr' , depending on its sign, must be carried 
away or supplied by phonons. 

There are three low-energy phonon modes in graphene: 
longitudinal, transverse, and flexural. In the long wave¬ 
length limit, the first two have the usual acoustic dis¬ 
persion, ojq = while the third one has a quadratic 
dispersion. For the same frequency, long-wavelength flex¬ 
ural phonon modes have a much larger wave-vector and, 
accordingly, reside in a much larger phase space. We 






4 



FIG. 3: UnpxiXjy) (a) and the ground state wave function (b) of a single hydrogen adatom above one of the central carbon 
sites of circumcoronene. The potential energy Uppp{x,y) was restricted to a disk of radius a/2 outside of which the energy 
was kept constant. The Shrodinger equation was then solved in a domain extending to the regions of the nearest neighbors. 
Only the region of the central carbon site is shown here. The contour lines in (a) are equally spaced by 0.1 eV while in (b) they 
are equally spaced by 0.5 nm“^. 


expect flexural modes to dominate phonon-assisted pro¬ 
cesses, similarly to how they dominate the phonon resis¬ 
tivity in electron transport!^. 

The flexural phonon’s quadratic spectrum follows from 
the Hamiltonian for the out-of-plane sheet oscillations, 
which, in the long-wavelength limit, i^^ 


Tfo = 


1 

2 


dh 


ph^ + D{\l‘^hf 


(5) 


with p and D denoting the two-dimensional mass den¬ 
sity and the flexural rigidity constant of graphene, re¬ 
spectively. The Hamiltonian ([^ describes unstrained 
graphene. In the presence of strain the Laplacian in Hq 
is replaced with Vh, in which case assumes the form 
of the usual membrane Hamiltonian. Below, we consider 
an unstrained sheet. As follows from the Hamiltonian 
(§, the spectrum of flexural modes is 

( 6 ) 

Standard quantization of the Hamiltonian ([^ leads to 
the following expression for the phonon flexural distur¬ 
bance 



energy of the atom’s oscillatory motion induced by the 
phonons. 


’Hn-ph = ( 8 ) 

Substitution of the phonon operator Q into this expres¬ 
sion yields the hydrogen-phonon coupling Hamiltonian 

Tdli—ph — Tdahs T Td em ^scf (9) 

which consists of three terms. The first term describes 
absorption of two phonons with wave-vectors q and k, 

nabs = - Y E (10) 

^ q.k 

where Aq = yJ‘hujq/2pA. The second term is the Hermi- 
tian conjugate to the first one, 

nem = ~Y. ( 11 ) 

^ q,k 

and describes the emission of two phonons. Finally, the 
term 


in terms of the phonon creation and annihilation opera¬ 
tors; A stands for the normalization area of the sheet. 

The hydrogen adatom represents a mass defect stuck 
to the graphene lattice. The coupling of the atom to the 
phonon held Q is described by the additional kinetic 


nsc = MY^ (12) 

q.k 

describes scattering processes, in which a phonon with 
wave-vector q is absorbed and a phonon with wave-vector 
k is emitted. 

In what follows we are going to utilize the perturba¬ 
tion theory to take into account the hydrogen-phonon 
coupling Such approach is justihed because of the 
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light mass of hydrogen M <C Me, which ensures that the 
coupling is going to be rather weak. 

The probability of phonon-assisted hydrogen hopping 
from site R to site R' follows from the Golden rule. For 
example, when hopping occurs with a decrease of the on¬ 
site energy, Tiujq = W-r, — Wr,' > 0, only phonon emission 
and scattering processes are allowed. The transition rate 
for the emission is 


2-k 


From the Hamiltonian (11), we obtain that the matrix 
element for the emission transition is 


|(R';q,k|Hen»|R)P = u:qUJk{l + Nq){l + Nk) 




(14) 


where Nq = (e^^’ — 1)“^ is the Bose-Einstein distribution 
with j3 = h/ksT. 

The total emission rate is found from Eqs. (13l-(14) 
by integrating over all phonon momenta. After simple 
calculation one finds. 


IF,, 


f, ^ Wo — W 

1287rp74 J 1 - e-/3- 1 - 
0 


Similarly, the probability of hydrogen atom hopping 
assisted by the absorption of a phonon with frequency 
oj and emission of a phonon with the higher frequency 
Wo -I- w is given by 


Wse 


-^^^R-R' f, ZU Wq-I-W 

32npD J e/5- - 1 ■l_e-/3(-o+-y- 
0 


The transition rates (15) and (16) can be easily calcu¬ 
lated in the limits of low and high temperatures. It turns 
out that the scattering channel is the dominant mecha¬ 
nism of phonon-assisted hopping when temperatures are 
high. On the other hand, emission prevails at low tem¬ 
peratures. This is expected since the number of phonons 
available for scattering is small in that case. 


facilitating tunneling than the emission transitions. The 
T^-dependence of the scattering-assisted probability can 
be understood as follows. Phonons involved in the tran¬ 
sition have energy of the order ksT. One power of tem¬ 
perature arises from the number of such phonons avail¬ 
able. The other two powers appear due to the fact that 
high-frequency phonons interact stronger with a hydro¬ 
gen adatom, since they result in a larger coupling Hamil¬ 
tonian 


B. Low temperatures, /Iwq S> 1 


In this limit, all exponentials with negative arguments 
can be discarded in both Eqs. (IT^ and (16). Addition¬ 


ally, since in the scattering channel the incident phonons 
have frequency w ^ 1//? <C wq, one may neglect w in the 
numerator of the integrand of Eq. (16). Using the fact 
that duu^/{e^ ~ 1) = we find, 


fFem 1 , . f WoV4, 

IF«, ]- 192ttpD '^°\{7rkBT/h)\ 


(18) 


The predominance of emission over scattering in the low- 
temperature limit is driven by the scarcity of phonons in 
the initial state. 


C. Hopping with the increase in energy, wq < 0 


When a hydrogen atom tunnels from a site with a lower 
interaction energy VFr to a site with a higher energy Wr> 
it needs to pick up the extra energy to do so. This is pos¬ 
sible by either absorbtion of two phonons via the process 
described by the Hamiltonian (10) or via the scattering 
processes already discussed above. The corresponding 
rates can be calculated from the same Golden rule for¬ 
malism. However, they also follow immediately from the 
detailed balance principle, which provides. 


Wabsi-UJo) = 

lF,,(-Wo) = lF,,(wo)e-/^‘"G (19) 


A. High temperatures, /3wo <C 1 


In this limit, the denominators in the emission proba¬ 
bility (15) can be expanded to the linear order in small 


(3. In the scattering probability ( |16[ ) one may set wq = 0 
since typical phonons participating in the process have 
much higher frequency, w ~ 1//3 ^ wq. Using the iden¬ 
tity, duv?/s\Rh^{u/2) = 47r^/3, we obtain, 


Wem \ fkBTV f wo/4, 

W,, j 32 tt pD \ n ) \^^kBT/3fi. 


As mentioned above, at high temperatures the 
scattering-assisted processes are much more efficient in 


At high temperatures, /3wo <C 1, therefore, the transi¬ 
tion rates for the upward (in energy) tunneling are vir¬ 
tually the same as for the downward transitions. At low 
temperatures, /3wo ^ 1, however, the upward transitions 
are strongly suppressed. 


IV. HYDROGEN-ELECTRON INTERACTION 

Inelastic tunneling of a hydrogen atom can also occur 
with the excess (or deficit of) energy transferred to (from) 
the conduction electrons. At low energies, the latter have 
the Dirac spectrum e = ±up. 

The amplitude of this process is proportional to the 
matrix element of the interaction of a hydrogen atom 
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with the electric charge of a conduction electron. A hy¬ 
drogen adatom may interact with conduction electrons 
via Coulomb forces. Since the adatom is neutral, this 
interaction, in the lowest order, results from the dipole 
moment associated with the carbon-hydrogen bond and 
has the form. 


ed|| ■ (r-R) 
k\y — RP 


( 20 ) 


where r is the coordinate of the electron and d|| is the 
adatom’s in-plane dipole moment. The coefficient k is the 
effective dielectric constant of graphene, which describes 
screening of static electric fields by other conduction elec¬ 
trons. 

From symmetry considerations, however, it is clear 
that the average dipole moment will be pointed in the 
direction perpendicular to the plane of graphene and thus 
the expectation value of the Hamiltonian ( |20[ ) in the 
ground state of a hydrogen adatom is zero, (0]ff'|0) = 0. 
As long as there is no spontaneous symmetry breaking, 
we need to explore the second order correction. 

This second order correction is analogous to the van 
der Waals interaction between two atoms where the in¬ 
teraction arises as a result of quantum fluctuations. The 
main difference comes from the fact that in our situation 
one of the particles (the electron) has a net charge, so 
that the resulting interaction falls of as l/r^, rather than 
as 1/r®. Indeed, the second order energy correction is 
finite and given by the standard expression. 


n 


Eq — En 


( 21 ) 


where the summation is taken over all excited states n. 

The actual energy levels and dipole moment’s matrix 
elements for a hydrogen atom sitting above a carbon 
atom in a graphene crystal can only be determined by 
first-principles calculations which are beyond the scope 
of the present paper. However, it is not difficult to esti¬ 
mate the correction (21) by the order of magnitude. The 
non-diagonal matrix elements are (n|d|[|0) ^ eas, where 
Os is the Bohr radius. Similarly, the difference Eq — 
is of the order of the Rydberg energy, e^ j2a b- As a re¬ 
sult, AA, which should be identified with the effective 
hydrogen-electron coupling energy, is. 


y(r-R) = AA = -C- 


2„3 


e a' 


Ir-R|4- 


( 22 ) 


where C is an unknown positive dimensionless constant. 


Expression (22) can also be obtained by purely clas¬ 


sical considerations. The electric field produced by the 
conduction electron, and acting on the hydrogen atom, is 
e/ftr^; for simplicity we set R = 0. This field is (r/os)^ 
times weaker than the typical atomic fields, It, 

therefore, leads to the displacement of the atomic elec¬ 
tron of the order, aB(aB/r)^, resulting in the induced 


dipole moment dy ^ ea^/r^, from which Eq. (22) fol¬ 
lows. Note that, for a free hydrogen atom, the dimen- 
sionless constant is known exactly^, C = 9/2, suggesting 
that C ~ 1 — 10 when the atom is sitting atop a graphene 
sheet. 

We can now write the hydrogen-electron interaction 
(22) in the second-quantized form with the help of the 
electron creation and annihilation operators. 


Hn-e = 2^ E E 

p,p' a,P 

X 4pC/3p'[l+ a/3cos(p,p')]. (23) 


In this expression the indices a and /3 assume two val¬ 
ues: -|-1 for electrons in the upper Dirac cone and —1 in 
the lower cone; the trigonometric factor in the brackets 
of the second line is coming from the pseudospin projec¬ 
tion of the initial electron state |/3p') onto the final state 

lap)- 

The Eourier transform of the interaction potential di¬ 
verges at low distances and should be cut-off there by the 
Bohr radius, 


V{h<i)= J d?re ^‘^'*'I/(r) » —2-kC ^ j ^Jo{qr) 


= —ttC 


e^as 


(24) 


The interaction is essentially short-range with the 
strength independent of the electron momenta. This is 
the result of the fast decay of the interaction ( [2^ with 
the distance and the fact that only long-wavelength elec¬ 
trons, pas <C can participate in the interactions, since 
energy involved in inelastic tunneling processes, Tiwo, is 
small compared with the electron bandwidth. 

The rate of inelastic tunneling transitions facilitated 
by the electron scattering is given by the Golden rule 
formula with the matrix element provided by Eq. (23): 


Wh-c = 


37rC2^4„2 


B t2 


ATfn 


5^2 


^R-R 


oo oo 

' E / “ ‘^ap )'n/3k 

0 0 

x6{huJo — avp + (dvk), (25) 


where Uap = [exp (avp/T) -\- 1]“^ is the Fermi-Dirac dis¬ 
tribution of electrons in the two cones. The overall co¬ 
efficient in Eq. ( |2^ takes into account the existence of 
two Dirac cones as well as the two-fold spin degeneracy. 
Note that [1 -I- a/3 cos (p, p')]^ yields the factor 3/2 upon 
the angle averaging. 

At high temperatures, ksT ^ Tiujq, the main contribu¬ 
tion comes from transitions that occur within the same 
cone, a = f3. In this limit huio can be neglected. After a 
simple integration we obtain. 


Wn-e 




(26) 
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At low temperatures, ksT hujQ, the dominant pro¬ 
cess results from the lifting of an electron from the lower 
cone, /3 = —1, into the upper cone, a = 1, yielding the 
following inelastic tunneling rate, 


WH-e = 


STi'^k'^v^ 


'R-R' 


(27) 


Comparing these results with the rates for phonon- 
assisted tunneling from Section we see that the elec¬ 
tron mechanism leads to the same temperature depen¬ 
dence as the dominant phonon processes: phonon scat¬ 
tering at high temperatures, see Eq. (171, and phonon 
emission at low temperatures, see Eq. (181. The reason 


for this coincidence is not difficult to understand. For 
example, at high temperatures electrons participating in 
the tunneling are confined in a smaller region of phase 
space; however, their coupling to hydrogen is stronger, 
since the coupling Hamiltonian (23) remains constant 


in the long-wavelength limit while coupling to phonons, 
Eq. (12), vanishes there. 


V. DISCUSSION 


Because of the similar temperature dependence, the 
ratio of the two rates takes a rather simple form. In¬ 
terestingly, the two mechanisms appear to lead to rates 
of the same order of magnitude. At high temperatures, 
ksT 3> Tiojq, the ratio of the two rates can be written as 


WH-e 


= 247r^C 


2^2 


(P^\ (^\ 


\?iKv) \ M ) \Mv'^) 


(28) 


The first ratio in this expression is /fiKV = 1.4, given 
the well-known value of the effective dielectric constant 
in graphen^^, k = 1 -|- jAfiv = 1.55, and the Fermi 
velocity, u = 1 x 10® m/s. The second ratio is of the 
same order, pa\/M = 16a^/(-\/3a^) = 1.3, and is easily 
calculated from the value of the lattice spacing, a = 0.14 
nm, and the fact that the graphene lattice is honeycomb. 
Given the gen erally quoted value for the flexural rigid¬ 
ity const aniPS D = 1.2 eV, the last ratio is very small, 
D/Mv^ = 1.4 X 10“^. However, this smallness is com¬ 
pensated by the large value of the numerical prefactor in 
Eq. (28). The constant C is not known. However, if we 


utilize its value for a free hydrogen atom in the ground 
stat^^, C = 4.5, we obtain that the ratio of the electron 
to the phonon contributions is 


Wn-e 

W^rr,^ 


1.45 


(29) 


At low temperatures, fcsT <C hujQ, the ratio is given 
by the expression similar to Eq. (28) where there is an 


additional coefficient 4, as follows from Eq. (27) and the 


first line of Eq. (18). As a result, we have Wn-elWem ~ 
5.83 there. 

The rates for inelastic hydrogen tunneling obtained 
in the preceding sections are rather small, predicting 


a very slow diffusion of hydrogen. Let us estimate 
the rate of hydrogen hopping in the dominant regime. 
Namely, consider the case of high temperatures where 
the scattering channel dominates, and the transition 
rate is given by the second line of Eq. 0- Using 
agaiiP^ D = 1.2 eV, and the overlap integral calcu¬ 
lated earlier, /r_r' ~ 3.5 x 10“®, we find that in the 
high temperature limit the total hopping rate provided 
by the electron-assisted and phonon scattering rates is 
W = Wh-b + Wsc ~ 1.2(T/To)^Hz in terms of the ref¬ 
erence temperature Tg = 300 K. 

It is also instructive to restate the hopping rate W in 
terms of the diffusion coefficient V. To estimate its order 
of magnitude we assume that ksT exceeds Tl\ujq \ on each 
site, so that the probability of tunneling W is uniform 
across the extent of the entire graphene sheet and is set 
by the temperature T. The diffusion coefficient for a 
particle hopping on a honeycomb lattice with a given 
probability W is then, 

X>=ilUa^ (30) 

At room temperature the diffusion coefficient is of the 
order V ~ 10“^^ m^/s. Given this value of the dif¬ 
fusion constant, it would take a hydrogen atom a time 
T ~ I'D ~ 100 years to diffuse across a micron-sized 
graphene flake. 

In this paper we disregarded polaronic effect. The lat¬ 
ter could be expected to renormalize mobility of a moving 
particle, by some numerical factor. On the other hand, 
the suppression of diffusion due to the effective interac¬ 
tion studied here appears to be a much stronger effect, 
by orders of magnitude. 

Because of the slow character of hydrogen adatom dif¬ 
fusion, hydrogen desorption and absorption might be¬ 
come important in specific situations. Those processes 
depend on the ambient conditions and are beyond the 
scope of this paper. 


VI. SUMMARY 


For electronics applications utilizing graphene samples 
doped with resonant adatoms, such as hydrogen, it is im¬ 
portant to predict how fast the adatom diffusion would 
occur at room temperatures. Fast diffusion could be fa¬ 
vorable for the device manufacturing, but it would be 
detrimental for the device longevity. On the other hand, 
fast diffusion would be required for the realization of var¬ 
ious phases resulting from adatom orderin^^^^^^. 

Our calculations show that the diffusion of hydrogen 
adatoms, the lightest dopants, is very slow. A single 
adatom placed on an ideal graphene sheet would propa¬ 
gate with the velocity of the order of Icm/s. However, 
the presence of other such adatoms suppresses adatom 
propagation dramatically. This is the result of a peculiar 
electron-mediated effective interaction between adatoms. 
















Two features make this interaction efficient in suppress¬ 
ing diffusion. First, it is a is long-range interaction, so 
that it is considerable even when adatoms are tens of 
nanometers away from each other. Second, it changes 
sign depending on whether adatoms reside on the same or 
opposite sublattices. In order to move, an adatom would 
have to hop to a nearest neighbor carbon site and change 
its energy in the process by as much as tens of meV. 
Such inelastic hopping requires assistance from phonons 
or electron-hole excitations. 

Let us emphasize that the suppression of the diffu¬ 
sion discussed in this paper is due to the graphene band 
structure featuring two Dirac points, which results in a 
quantum interference of the electron band states prop¬ 
agating over two sublattices in the honeycomb arrange¬ 
ment of carbon atoms. The effective hydrogen-hydrogen 
interaction, which has opposite sign on the two sublat¬ 
tices, originates from coupling of hydrogen adatoms to 


such band electrons. 

As a result, we find a single adatom hopping time to 
be in the range of milliseconds at room temperature. 
This indicates that the diffusion is sufficiently slow that 
hydrogen-doped graphene devices are feasible. However, 
it is too slow for any gate-control of the hydrogen adatom 
distribution to be useable in fast switching devices. 
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